#include <stdio.h>
#include <gmp.h>
#include <time.h>
#include <stdlib.h>
#include <pthread.h>

#define FOURBYTES 4294967296

mpz_t Gmax_it, Gcount, Git;                           //numero total de iteracoes, soma de todos os pontos dentro do circulo e numero de iteracoes por thread (nIteracoes/nthreads)
pthread_mutex_t mutex;

void *work(){                                         //funcao a ser passada para cada thread trabalhadora
	mpz_t it, count;                              //contador local iteracao e de pontos dentro do circulo (individual em cada thread)
	mpz_inits(it, count, NULL);                   
	gmp_randstate_t state;                        //variavel que guarda o estato do gerador de numeros aleatorios
	double x, y;                                  //coordenadas x,y do ponto sorteado

	gmp_randinit_mt(state);                       //inicializacao e semente do gerador de numeros aleatorios
	gmp_randseed_ui(state, time(NULL));

	for(mpz_set_ui(it, 0); mpz_cmp(it, Git); mpz_add_ui(it, it, 1)){	//for(it = 0; it != Gmax_it; it++);
	
		x = (double) (FOURBYTES - gmp_urandomb_ui(state, 32)) / (double) FOURBYTES; //sorteia x e y
		y = (double) (FOURBYTES - gmp_urandomb_ui(state, 32)) / (double) FOURBYTES; 	
		if( x*x + y*y <= 1.0 )                                                      //verifica se P(x,y) esta dentro do circulo de raio 1
			mpz_add_ui(count, count, 1);  //count++
	}

	pthread_mutex_lock(&mutex);		      //regiao critica: cada thread soma nas variaveis globais os resultados que obtiveram
		mpz_add(Gcount, Gcount, count);       //numero de pontos dentro do circulo
	pthread_mutex_unlock(&mutex);		      //fim da regiao critica

	mpz_clears(it, count, NULL);                  //limpa variaveis locais
	gmp_randclear(state);                         //limpa estado do gerador de numeros aleatorios local
	pthread_exit(NULL);                           //fim da thread
}

int main (int argc, char *argv[]){
	
	pthread_t *threads;                           //threads
	unsigned int nthreads, i;
	int erro;             			      
	FILE *file = fopen("PImcarloParalelo", "w+");
	mpf_t PI, mpf_count, mpf_it;                  //variavel para guardar o PI e floats temporarios para efetuar a divisao
	
	mpf_set_default_prec(35000000);               //inicializa variaveis com precisao de 35 milhoes de bits
	mpf_inits(PI, mpf_count, mpf_it, NULL);
	mpz_inits(Gmax_it, Gcount, Git, NULL);

	if(argc != 3){
		printf("Passe por linha de comando o numero de iteracoes seguido do numero de threads\n");
		return 0;
	}

	gmp_sscanf(argv[1], "%Zd", Gmax_it);          //pega argumentos da linha de comando
	sscanf(argv[2], "%d", &nthreads);

	threads = (pthread_t*) malloc (nthreads * sizeof(pthread_t));  //o usuario pode definir o numero de threads, que sao alocadas dinamicamente
	pthread_mutex_init(&mutex, NULL);	                       //inicializacao do mutex
	mpz_tdiv_q_ui(Git, Gmax_it, nthreads);                     //divisao do numero de iteracoes pelo numero de threads
	
	if(file == NULL){
		printf("Erro ao abrir o arquivo\n");
		return 0;
	}
	
	printf("Iniciando calculo do pi com %s iteracoes pelo metodo de Monte Carlo paralelo usando %s threads.\n", argv[1], argv[2]);

	for(i = 0; i < nthreads; i++){
		erro = pthread_create(&threads[i], NULL, work, NULL);  //cria threads e inicializa retornando um codigo de erro se houver
		if(erro){ 
			printf("Erro ao criar thread %d \n", i);
			exit(-1);
		}
	}
	
	for(i = 0; i <nthreads; i++){
		erro = pthread_join(threads[i], NULL);                 //aguarda a finalizacao de todas as threads, retornando um codigo de erro se houver
		if(erro){
			printf("pthread_join retornou erro: %d\n", erro);
			exit(-1);
		}
	}

	mpf_set_z(mpf_count, Gcount);                                  //conversao dos contadores para float, para que se possa fazer a divisao
	mpf_set_z(mpf_it, Gmax_it);
	
	mpf_div(PI, mpf_count, mpf_it);                                //PI = (double) (4*count) /(double) it;
	mpf_mul_ui(PI, PI, 4);                                         
	
	printf("Calculo finalizado.\nEscrevendo resultado em PImcarloParalelo\n");	
	
	mpf_out_str(file, 10, 0, PI);                                  //grava resultado no arquivo, em decimal com precisao padrao
	
	fclose(file);
	file = NULL;

	mpf_clears(mpf_count, mpf_it, PI, NULL);                       //libera variaveis
	mpz_clears(Git, Gmax_it, Gcount, NULL);
	free(threads);	
	pthread_mutex_destroy(&mutex);
	pthread_exit(NULL);

return 0;
}


